Any Progress on inverse pyrolysis modelling

with ensemble learning methods?

Patrick Lauer
University of Wuppertal

Content

What is pyrolysis?

Burning of solids can be separated in two phases:

  • Thermochemical decomposition of solid material and phase change from solid to gas phase (Pyrolysis)
  • Chemical reaction in the gas phase (Combustion)

To predict fire spread, we need to model burning of solids, hence pyrolysis.

Pyrolysis

pyro

Figure 1. [Schematic of pyrolysis [10]]

How do we model pyrolysis?

Boundary condition

(1)
\[-k_{s,1} \frac{T_{s,1}^{n+1}-T_{s,0}^{n+1}}{\delta x_{\frac{1}{2}}}=\dot{q}_c'' + \dot{q}_r'' \]

Heat conduction

(2)
\[\rho_s c_s \frac{\partial T_s}{\partial t} = \frac{\partial}{\partial x} (k_s \frac{\partial T_s}{\partial x})+\dot{q}_s''' \]

Reaction rate:

(3)
\[r = A Y^n \cdot e^{-\frac{E_a}{RT}} \]

Parameter overview

Parameter
Activation energie $E_a$
Pre-exponential factor $A$
Reaction order $n$
Density $\rho$
Conduction coefficient $k$
Heat capacity $c$
Emissivity $\epsilon$
Heat of reaction $\Delta H$

How do we get these parameters?

Find parameters with small scale experiments and mathematical fitting, scale up to parts and devices

Usual experiments:

  • Thermogravimetrical analysis
  • Cone calorimeter
  • Micro combustion calorimeter
  • …

Exkurs: Thermogravimetrical analysis (TGA)

A small specimen (mg scale) gets heated in a furnace with a constant or transient heating rate. Mass loss of the specimen is captured. It allows to estimate reaction kinetics of this sample.

tga

Figure 2. TGA

Approaches

  • Forwad fitting
    • Basic graphical fitting [4, 5, 8]
    • Advanced automated fitting [3]
  • Inverse modeling [9]
    • Optimization algorithms [1, 7]
    • Machine learning

Forests of randomized trees

Common Methods:

  • Random Forests [2]
  • Extremly Randomized Trees [6]

Concept

  • Ensemble learning methods

Advantages

  • Efficient for big data sets
  • Fast to train
  • Easy to parallelize

Method

  • Train a model to predict reaction kinetic parameters with given reaction rate
  • Case study: mockup TGA experiment with constant heating rate
  • All data used is randomly generated with the pyrolysis model

Method II

invers

Figure 3. Invers modelling then

Method III

invers2

Figure 4. Invers modelling now

Process I

  1. Generating sample data set with the regarding model
    • Samples for 1, 2 and 3 reactions, with 3 heating rates applied each
    • Up to 1M samples generated, with $r(T)$ and $A_n$, $E_n$

Exkurs: Sampling I

Problem:

  • Arrhenius equation is an exponential function
  • $A_n$ might be $10^{10}...10^{40}$
  • $E_n$ might be $10^{1}...10^{20}$
  • Hence we can't sample from a uniform distribution

Challenge:

  • Find a distribution to sample from

Exkurs: Sampling II

Solution:

  • Introducing $T_p$ (reference temperature), $\Delta T_p$ (reference range) and $r_p$ (reference rate), characterizing a triangle of temperature at maximum mass loss rate and width of the peak

peak

Figure 5. $T_p$, $\Delta T_p$, $r_p$

Exkurs: Sampling III

  • $T_p$ and $r_p$ are sampled from a uniform distribution and then mapped to $A_n$ and $E_n$
  • Mapping between $A_n$, $E_n$ and $T_p$, $\Delta T_p$, $r_p$ is done with these equations:
(4)
\[E_{i,1} = {\frac {er_{p,i}}{Y_{s,i}(0)}}{\frac{R{T}_{p,i}^2}{\dot T}} \]
(5)
\[A_{i,1} = {\frac {er_{p,i}}{Y_{s,i}(0)}} e^{\frac {E}{RT_{p,i}}} \]
(6)
\[\frac {r_{p,i}}{Y_{s,i}(0)}={\frac{2 \dot{T}}{\Delta T}}(1-\nu_{s,i}) \]

Process II

  1. Splitting data set in two independent sets (75 % training data set and 25 % validation data set)
  2. Train model with training data set
    • Input: $r(T)_{train}$
    • Output: $A_{n, train}$, $E_{n, train}$
    • Model adopts to transform input to output
  3. Validate trained model by feeding $r(T)_{prescribed}$ of validation data set and check for expected outcome

Process III

  1. Recalculate $r(T)_{predicted}$ with $A_n$, $E_n$, calculate RMSE between $r(T)_{validation}$ and $r(T)_{predicted}$

Process IV

  1. Evaluate
  2. Repeat with different algorithms and different hyperparameter settings

Results

for inverse replacement models with

  • 1 reaction, 3 heating rates,
  • 2 reactions, 3 heating rates and
  • 3 reactions, 3 heating rates

built with Extremly Randomized Trees algorithm

1 Reaction, 3 Heating rates

1r

Figure 6. Histogram of RMSE Test vs. Predicted Data (12.5k cases)

2 Reactions, 3 Heating rates

2r

Figure 7. Histogram of RMSE Test vs. Predicted Data (250k cases)

2 Reactions, 3 Heating rates II

example

Figure 8. Random example with ET

3 Reactions, 3 Heating rates

3r

Figure 9. Histogram of RMSE Test vs. Predicted Data (250k cases)

Evaluation

Advantages over other methods

  • Trained model is fast (instant result)
  • Trained model is portabel
  • Results are pretty good (by now in 53 % of presented case)
    • better for 1 reaction
    • worse for 3 reactions
  • If no perfect fit was found, it is at least a good starting point for other methods

Disadvantages over other methods

  • Generating samples is costly
  • Training a model is costly
  • Results are only good in 53 % of presented case

Performance assesment against evelutionary optimization algorithms

SCEUA ET
Setup time - 15…150 min
Calculation time 15…150 min seconds

Once an inverse replacement model with ET is set up, it is advantageous over SCEUA, since it gives instant results.

Performance assesment against training data set I

Question:

  • Is the training data set itself large enough to contain a good fit for the data that shall be predicted?

Process:

  • Split data set in training and testing (75/25)
  • Calculate RMSE for each training case with each testing case
  • Save best fit value for each testing case and plot
  • Drawback: High computational costs, since $1.87 \cdot 10^{11}$ RMSE calculations have to be conducted for a dataset with $10^{6}$ experimental sets
  • So in this case just 10 % of the testing data set is used and compared to the training data set: only $1.87 \cdot 10^{10}$ calculations

Performance assesment against training data set II

eval

Figure 10. Histogram of RMSE Test vs. Training Data (25k cases)

Performance assesment against training data set III

eval

Figure 11. Histogram of RMSE Test vs. Training Data

eval2

Figure 12. Histogram of RMSE Test vs. Predicted Data (250k cases)

Discussion I:

Ugh! What happend here?

Possible explanations:

  • Some error in the code?
  • My solution actually performs that bad?
  • Ill posed problem and “wrong” solutions found from database?

Discussion II:

surface

Figure 13. Sensitivity of log(A) and E changes to RMSE

Discussion: Evaluation

  • When calculating the fit with $RMSE$, we get some numbers that represent the fit
  • Without normalization it's not even comparable to other cases
  • Even with normalization, question still is:
    • What is a good fit?
      • either as a numeric value
      • or utilizing other criteria to assess
  • Let's discuss!

Next steps I

  • Use larger sample data sets
    • implementing a database (HDF5)
  • Compare to other machine learning models

Next steps II

  • Couple with Heat conduction
  • Try with different Pyrolysis models
    • gpyro
  • Validate with real data
    • first step: PMMA from MacFP

Progress on PROPTI

invers

Figure 14. Inverse Modelling

PROPTI workflow

workflow

Figure 15. PROPTI workflow

PROPTI

qr

https://github.com/FireDynamics/propti

New functionalities

Focus: Cost Function

  • determines the deviation between target data and a model response, e.g. between experimental data and simulation results.

Several different cost functions are presented for estimating material parameter sets that allow the simulation of pyrolysation of solid polymers.

Cost Function - Single Point

sp

$\text{NRSE}=\frac{|\hat{y}_t-y_t|}{y_t}$

Cost Function - Threshold

th

$\text{THR}_{min}= \min\Big\{|\{t | \hat{y}(t) > y(t_0)\}-t_0|, \max\{|t_{min}-t_0|, |t_{max}-t_0|\}\Big\}/t_0$

Cost Function - RMSE

rmse

$\text{RMSE}= \frac{1}{y_N \sqrt{n_y}} \sqrt{\sum_{t=1}^{n_y}(\Delta y_t)^2} \text{ with } \Delta y =\hat{y}_t-y_t$

Cost Function - RMSE BANDS

rmseband

$\text{RMSE}= \frac{1}{y_N \sqrt{T}} \sqrt{\sum_{t=1}^T(\Delta y_t)^2} \text{ with } \Delta y = \begin{cases} 0 & y_{t, lb} \leq \hat{y} \leq y_{t, ub}\\ \hat{y}_t-y_{t, lb} & \hat{y}_t < y_{t, lb} \\ \hat{y}_t-y_{t, ub} & \hat{y}_t > y_{t, ub} \\ \end{cases}$

Cost Function - RMSE RANGE

rmserange

$\text{RMSE}= \frac{1}{y_N \sqrt{T}} \sqrt{\sum_{t=1}^T(\Delta y_t)^2} \text{ with } \Delta y = \begin{cases} 0 & (1-r)y_{t} \leq \hat{y} \leq (1+r)y_{t}\\ \hat{y}_t-(1-r)y_{t} & \hat{y}_t < (1-r)y_{t} \\ \hat{y}_t-(1+r)y_{t} & \hat{y}_t > (1+r)y_{t} \\ \end{cases}$

Cost Function - Combination

combined

$E=\sum_{i=0}^I (w_i \cdot \text{RMSE}_i) + \sum_{j=0}^J (w_j \cdot \text{THR}_j) + \sum_{k=0}^K (w_k \cdot \text{NRSE}_k)$

Discussion

A cost function that uses an area as a target, provides means to incorporate the uncertainty observed in the experiments.

RMSE requires exact matches of the data points, while slight variations in the other cases could still fall inside the target area. BANDS and RANGE could be useful to account for variance that is encountered when repeating a single experiment multiple times and allow for its representation during the IMP.

The ability to combine cost functions in different ways allows to target multiple values, like heat release rate or surface temperature, and their unique features, like heat release peaks on different experimental setups (e.g. different heat fluxes or gas atmospheres), as these may be of crucial importance for the real-scale applications, especially for flame spread modelling.

References

[1]Lukas Arnold, Tristan Hehnen, Patrick Lauer, Corinna Trettin, and Ashish Vinayak. “PROPTI–A Generalised Inverse Modelling Framework.” In Journal of Physics: Conference Series, 1107:032016. IOP Publishing. 2018. 🔎
[2]Leo Breiman. “Random Forests.” Machine Learning 45 (1). Springer: 5–32. 2001. 🔎
[3]Morgan C Bruns, and Isaac T Leventon. “Automated Fitting of Thermogravimetric Analysis Data.” In . Interflam. 2019. 🔎
[4]Joseph H Flynn, and Leo A Wall. “A Quick, Direct Method for the Determination of Activation Energy from Thermogravimetric Data.” Journal of Polymer Science Part C: Polymer Letters 4 (5). Wiley Online Library: 323–328. 1966. 🔎
[5]Henry L Friedman. “New Methods for Evaluating Kinetic Parameters from Thermal Analysis Data.” Journal of Polymer Science Part C: Polymer Letters 7 (1). Wiley Online Library: 41–46. 1969. 🔎
[6]Pierre Geurts, Damien Ernst, and Louis Wehenkel. “Extremely Randomized Trees.” Machine Learning 63 (1). Springer: 3–42. 2006. 🔎
[7]Mohammed M. Khan, Archibald Tewarson, and Marcos Chaos. “Combustion Characteristics of Materials and Generation of Fire Products.” In SFPE Handbook of Fire Protection Engineering, edited by Morgan J. Hurley, Daniel T. Gottuk, John R. Hall Jr., Kazunori Harada, Erica D. Kuligowski, Milosh Puchovsky, Jose´ L. Torero, John M. Watts Jr., and CHRISTOPHER J. WIECZOREK, 1143–1232. Springer New York, New York, NY. 2016. doi:10.1007/978-1-4939-2565-0_36. 🔎
[8]Anna Matala, Chris Lautenberger, and Simo Hostikka. “Generalized Direct Method for Pyrolysis Kinetic Parameter Estimation and Comparison to Existing Methods.” Journal of Fire Sciences 30 (4). Sage Publications Sage UK: London, England: 339–356. 2012. 🔎
[9]Tatenda Nyazika, Maude Jimenez, Fabienne Samyn, and Serge Bourbigot. “Pyrolysis Modeling, Sensitivity Analysis, and Optimization Techniques for Combustible Materials: A Review.” Journal of Fire Sciences 37 (4-6). SAGE Publications Sage UK: London, England: 377–433. 2019. 🔎
[10]Guillermo Rein. “From Pyrolysis Kinetics to Models of Condensed-Phase Burning.” In . 2008. 🔎